sigma = c(0.2, 0.25, 0.18)
rho = matrix(c(1, 0.8, 0.4, 0.8, 1, 0.3, 0.4, 0.3, 1), byrow = TRUE, nrow = 3)

sigma2 <- (sigma %*% t(sigma)) * rho


w <- c(-0.5, -0.5, 1)
mu = c(0.1, 0.15, 0.12)
(mu.P <- t(w) %*% mu)

set.seed(1)
Data <- data.frame(Port1 = rnorm(1000, -0.5, 0.2),
                   Port2 = rnorm(1000, -0.5, 0.25),
                   Port3 = rnorm(1000, 1, 0.18))


(mean.R <- apply(Data, 2, mean))
(cov.R <- cov(Data))
(sd.R <- sqrt(diag(cov.R)))

(Sigma.P <- (t(w) %*% sigma2 %*% w))^0.5

library(quadprog)

Amat <- cbind(rep(1, 3), mean.R)  ## set the equality constraints matrix
mu.P <- seq(min(mean.R - 5e-04), max(mean.R + 
                                       5e-04), length = 300)  ## set of 300 possible target portfolio returns
sigma.P <- mu.P  ## set up storage for std dev's of portfolio returns
weights <- matrix(0, nrow = 300, ncol = ncol(Data))  ## storage for portfolio weights
colnames(weights) <- colnames(Data)



for (i in 1:length(mu.P)) {
  bvec = c(1, mu.P[i])  ## constraint vector
  result = solve.QP(Dmat = 2 * cov.R, 
                    dvec = rep(0, 3), Amat = Amat, 
                    bvec = bvec, meq = 2)
  sigma.P[i] = sqrt(result$value)
  weights[i, ] = result$solution
}

par(mfrow = c(1, 1))
plot(sigma.P, mu.P, type = "l", xlim = c(0, 
                                         max(sd.R) * 1.1), ylim = c(0, max(mean.R) * 
                                                                      1.1), lty = 3, lwd = 3)  ##  plot
## the efficient frontier (and
## inefficient portfolios below the
## min var portfolio)
mu.free = 1.3/253  ## input value of risk-free interest rate
points(0, mu.free, cex = 1, pch = "+")  ## show risk-free asset



sharpe = (mu.P - mu.free)/sigma.P  ## compute Sharpe's ratios
ind = (sharpe == max(sharpe))  ## Find maximum Sharpe's ratio
options(digits = 3)
lines(c(0, 2), mu.free + c(0, 2) * (mu.P[ind] - 
                                      mu.free)/sigma.P[ind], lwd = 4, lty = 1, 
      col = "blue")
## show line of optimal portfolios
points(sigma.P[ind], mu.P[ind], cex = 4, 
       pch = "*")  ## show tangency portfolio
ind2 = (sigma.P == min(sigma.P))  ## find the minimum variance portfolio
points(sigma.P[ind2], mu.P[ind2], cex = 2, 
       pch = "+")  ## show min var portfolio
ind3 = (mu.P > mu.P[ind2])  ## finally the efficient frontier
lines(sigma.P[ind3], mu.P[ind3], type = "l", 
      xlim = c(0, max(sd.R) * 1.1), ylim = c(min(mean.R) * 
                                               1.05, max(mean.R) * 1.1), lwd = 3, 
      col = "red")  ##  plot the efficient frontier
text(sd.R[1], mean.R[1], "Port1", cex = 1.15)
text(sd.R[2], mean.R[2], "Port2", cex = 1.15)
text(sd.R[3], mean.R[3], "Port3", 
     cex = 1.15)